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Abstract 

Although the existence of cosmic rays with energies extending well above 10 19 eV has 
been confirmed, their origin remains one of the most important questions in astro- 
particle physics today. Several source models have been proposed for the observed 
set of Ultra High Energy Cosmic Rays (UHECRs). Yet none of these models have 
been conclusively identified as corresponding with all of the available data. One 
possible way of achieving a global test of anisotropy is through the measurement of 
the information dimension, D\, of the arrival directions of a sample of events. D\ 
is a measure of the intrinsic heterogeneity of a data sample. We will show how this 
method can be used to take into account the extreme asymmetric angular resolution 
and the highly irregular aperture of a monocular air- fluorescence detector. We will 
then use a simulated, isotropic event sample to show how this method can be used 
to place upper limits on any number of source models with no statistical penalty. 

Key words: cosmic rays, anisotropy, fractal dimensionality, Centaurus A, dark 

matter halo, air fluorescence 
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1 Introduction 

The observation of Ultra High Energy Cosmic Rays (UHECRs) has now 
spanned nearly three decades. Over that period, many different source models 
have been proposed to explain the origin of these remarkable events. Recently, 
the Akeno Giant Air Shower Array (AGASA) reported clustering at small 
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angular scales for the events that were observed above 4 x 10 19 eV [1]. How- 
ever, this result could not be confirmed by the High Resolution Fly's Eye 
(HiRes) air fluorescence detector [2] despite the fact that Hires-I's monocu- 
lar exposure was more than twice that of AGASA within the pertinent energy 
range [3] . AGASA has also reported a possible excess of events in the direction 
of the galactic center for events with energies around 10 18 eV [4]. However, 
HiRes reported that it did not see anisotropies when examining harmonics 
in right ascension, a priori determined point sources, or enhancement in the 
supergalactic plane. Furthermore, prior analysis by the the original Fly's Eye 
showed no evidence of anisotropies when dependencies in galactic latitude 
and longitude, harmonics in right ascension, excess maps, and specific a priori 
determined point sources were examined [5,6]. However, Fly's Eye did find 
a small apparent excess in the direction of the galactic plane. [7]. In 1995, 
it was reported by Stanev et al. that the combined data of Haverah Park, 
Yakutsk, and AGASA showed an excess along the supergalactic plane with 
several potential point sources for events above 2 x 10 19 eV [8]. 

These conflicting reports call for developing a more global way in which 
one could determine if a given sample possesses any statistically significant 
anisotropy. We will show that by considering the information dimension, Dj, 
of a given sample, one can simultaneously look for anisotropies at all angular 
scales greater than the angular resolution of the sample by considering the 
intrinsic heterogeneity of that particular data sample. This method is quite 
robust in that it can easily accommodate both asymmetric angular resolutions 
and irregular apertures. Furthermore, in the event that a sample is shown to 
be consistent with an isotropic distribution, this method can be used to place 
upper limits on possible source models. Because only a single measurement 
is taken of the actual data, any number of source models can be considered 
simultaneously without incurring statistical penalties. 



2 Calculating Fractal Dimensions of a Data Sample 

Fractal dimensionality is a simple measure of the scaling symmetry of a struc- 
ture. By measuring the fractal dimension of a data sample, one can examine 
its self-consistency at different levels of magnification. There are several ways 
of exploiting this idea. From a computational perspective, the simplest is to 
use box-counting. For the most general case, the capacity dimension, D c [9], 
one partitions the sample space into equi-sized and equi-shaped "boxes" with 
edge size e: 




(1) 
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Here, N(e) is the minimum number of "boxes" with edge size e necessary to 
cover one's sample. 



The capacity dimension has a serious limitation: It only looks for the presence 
of the sample within the available space and does not consider variations in 
the density of the sample at a given point in space. In cases where the density 
may differ within the sample space, the appropriate alternative is to use the 
information dimension, D\ [10,11]: 

, = -£ Um ^f°sf-M (2 ) 



where -Pi(e) is the probability of finding a data point in the i-th box of edge 
size e. This is a particularly suitable measurement when considering a data 
set consisting of UHECR arrival directions with finite angular resolution. 

It should be noted that Dq and D\ are both particular cases of the g-dimension 
[12,13], 

D q = — - — lim (3) 

q 1-q e^o+ log 1/e ' y ; 



where 



N 



i=l 



We can then see that Dq = lim^o D q and that D\ = lim^i D q . 



3 Application to Arrival Direction Distributions for UHECRs 



In principle, it is simple to calculate the information dimension for a given 
sample of events. However, two complications arise when considering a set of 
arrival directions of UHECRs. First, the event directions are not known with 
complete precision. This makes the determination of Dj as e — > meaningless. 
Secondly, the determination of Di requires that the sample space be divided 
into equi-sized and equi-shaped bins. For a spherical surface, this is simply not 
possible. Nevertheless, there are workable solutions for both of these problems. 



3 




We will consider a hypothetical monocular air-fluorescence detector. This de- 
tector is largely based upon the characteristics of the Hires-I detector [14,15]. 
Our detector observes events with an angular resolution that is described by 
a highly asymmetric 2-d Gaussian. For a monocular air fluorescence detec- 
tor, angular resolution consists of two components: a±, in the determination 
of the angle, ip, within the plane of reconstruction, and 02, in the estimation 
of the plane of reconstruction itself. Figure 1 illustrates how this geometry 
would appear with a particular plane of reconstruction and a particular value 
for if). Intuitively, we can see that we should be able to determine the plane 
of reconstruction quite accurately. However, the value of ip is more difficult 
to determine because it is dependent on the precise results of the monocular 
reconstruction [14,15]. 

The actual parameterizations of o\ and <t 2 assumed are as follows: 

01 = 20°e" 1 - 51og i» £EcV + 4° (5) 



and 



a 2 = 100°e-°- 5A * + 0.4°. 



(6) 



Here, i?EeV is the primary energy of the shower in EeV. For the purpose of this 
study, the energy will be allowed to vary between 10 18 - 5 eV and 10 20 eV with 
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a differential spectral index of —2.7. In this scenario, a shower with a primary 
energy of 10 18 ' 5 eV will have U\ = 13.4°, while a shower with a primary energy 
of 10 20 eV will have <j\ = 5.0°. This difference can be attributed to the fact 
that larger showers have better defined profiles and a better signal-to-noise 
ratio. 

The factor, Ax, in equation 6 is the angular track length (in degrees) of the 
shower as observed by the detector, which is allowed to vary between 8° and 
30°. A shower with an observed track length of 8° will have a 2 = 2.2°, while a 
shower with an observed track length of 30° will have a 2 = 0.4°; a longer track- 
length leads to a more accurate determination of the plane of reconstruction. 
The distribution of o 2 values is largely independent of energy because while 
higher energy showers do lead to more longitudinal development, they are 
also brighter, which allows one to observe them at greater distances. These 
competing factors lead to A\ distributions that are virtually identical across 
the observed spectrum. 

In general, it should be noted the distributions of <i\ and a 2 values are relatively 
insensitive to the differential spectral index that is chosen. We ascertained this 
by considering two sets, one with a differential spectral index of —2.5 and one 
with a differential spectral index of —3.5. Even for a variation that was much 
larger than the accepted range of experimental values for the the UHECR 
spectrum [14-17], the value of &1 increased by only 11%. The value of <j 2 
remained unchanged. This can be explained by realizing for a steeply falling 
spectrum, the overwhelming majority of observed events in either case will 
occur in the first half decade of the measurement. This is a very small effect 
compared to the expected statistical fluctuations that would occur between 
two consecutive sets of observations. 

For the purpose of calculating Di, we can treat the arrival direction of each 
individual shower as a two dimensional elliptical Gaussian distribution with 
the parameters 01, cr 2 . The size of a bin's edge, e, will be allowed to take on a 
series of values, A9, which will be in an interval corresponding to a scale length 
of the sample or in the case of a smooth distribution, the smallest value that is 
computationally feasible. For finite events samples, we will use A6 ~ 0.5°. In 
the case of smooth distributions we will use a computationally limited value 
of Ad = 1/6°. The number of points in each shower direction distribution, 
-^Dist, will be determined by the mean value, <n ; >, necessary to assure that 
the fractional Gaussian fluctuations of the count, n;, in each bin, do not on 
average, exceed a predetermined value (i.e. for 5% fluctuations, ~ 500). 
For each value of e, the probability, -Pi(e), for the i-th bin will be: 



We calculate D\ for each value of e: 



We then determine D\ to be <Dj(e)> over the specified interval of e values. 
3.2 Latitudinal Binning 



For the purpose of calculating Dj, it is necessary that all bins be equi-sized 
and equi-shaped as we vary the size of the side of the bins, e. While it is 
impossible to achieve completely this criterion on the surface of a sphere, we 
will be able to do so approximately by adopting a latitudinal binning scheme. 

Latitudinal binning is achieved by first dividing the sky into N$ declinational 
(5) bands where each band has a width 



For each declinational band, the sky is then divided into N a> s bins in right 
ascension (a) where: 



N, 



a,8 



2-7T cos S dS 



2{N S ) 2 cos SdS 



7Y 



(10) 



The solid angle, AQ$ of each bin (in steradians) is: 



An, = *L^5!*f», ( „: 



with a minimum value of (A^) 2 (at the equator) and a maximum value of 
|(A#) 2 (at the poles) regardless of the value of N$. This provides us with 
bins that are all almost the same area and nearly square-shaped (with the 
exception of three triangular bins at each pole). The total number of bins in 
the sky can be approximated by: 

N sky ^4n(^) 2 = -(N 5 ) 2 . (12) 

71 71 



In figure 2, we visualize the latitudinal binning technique for a series of differ- 
ent N$ values. 
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Fig. 2. Hammer- Aitoff projection of latitudinal bins for four different values of 
N 5 — (a): N 5 = 5; (b): N s = 12; (c): N s = 30; (d): N 5 = 90. 

3.3 Application to the Calculation of Di 



We can now apply the preceding machinery to the calculation of D\: First, we 
need to normalize the event count in each bin by its respective bin area, AQ$: 



A^Dist -^Shower A Qfi ^Dist -^Shower (Ns) 2 A£ls 



(13) 



If we then realize that e = j^, we can obtain 



N 



^P^logPiiNs). 



(14) 



This expression is reminiscent of the the general formula for entropy from 
statistical mechanics: 



S = -kJ2Pr\ogp r ; 



(15) 



where p r is the probability of a particle being the r-th state and k is the 
Boltzmann constant, which can be thought of as a scaling constant based 
upon the intrinsic scale of the given particle. The information dimension, D\, 
is an analogous measurement of the heterogeneity of a given data set. 
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4 Calculating Di for Exposures of Different Source Models 

4-1 Exposure-Independent Source Descriptions 

We began by examining four different source models independently of detector 
exposure. This allows us to calculate the value of the information dimension, 
Di, without consideration to the detector aperture or statistical fluctuations 
from a finite event sample. The source models are: an isotropic source model, a 
dipole source model, a model with seven sources superimposed on an isotropic 
background, and a dark matter halo source model. 

4-1.1 Isotropic Model 

The first model that we will consider is an isotropic source model with distri- 
bution: 

^isotropic 

4.1.2 Dipole Model 

The second is the Centaurus A dipole source model first proposed by Farrar 
and Piran [18]. This model has a distribution of arrival directions characterized 
by a scaling parameter, a, which can take on any value between —1 and +1 and 
by 9, which is the opening angle between a given event arrival direction and 
the center of the dipole distribution at Centaurus A. The overall distribution 
is: 

^dipoic = 1 + acos6>. (17) 
4-1.3 Discrete Source Model 

The third model that we will consider is one with seven discrete sources su- 
perimposed on an isotropic background. For simplicity's sake, we will assume 
that all seven sources have an equal intensity indirectly determined by a pa- 
rameter, F s , which will be defined as the fraction of the entire event sample 
which originates in the seven sources. We will define our source direction to be 
the centroids of the seven hypothetical point sources proposed by the AGASA 
collaboration [1]. The equatorial coordinates used for each point source are 
listed in table 4.1.3 

The arrival directions for the events from each source are assumed to be sub- 
jected to magnetic smearing. That is, in the course of traveling through space 
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Table 1 

Coordinates used for the centers of seven discrete sources. These coordinates cor- 
respond to the centers of the seven clusters reported by the AGASA Collaboration 
[!]■ 

from the source to the point of observation, the velocity vector of the event 
is subject to bending in the galactic and extra-galactic magnetic fields. We 
assumed that this bending produces an apparent source that can be charac- 
terized by a Gaussian distribution: 

A9 (A6>) 2 

P(A9) = ^e^, (18) 



where P(A9) is the probability that an event will be observed with an opening 
angle, A9, from the nominal direction of the apparent source. We will also 
assume that the arrival directions of the events from all the sources are subject 
to the same degree of magnetic smearing, parameterized by A = 1.5105 • A. 
The parameter, A , corresponds to the 68% confindence interval in A9. For 
this paper, we will assume that A = 5°. It should be noted that the nominal 
direction of the apparent source is not necessarily the direction of the actual 
source because the possibility exists that the path of the events in question 
traveled through large regions of homogeneous magnetic fields. 



4.1.4 Dark Matter Halo Model 

The fourth model that we will consider is a dark matter halo source model. 
Dark matter halos are characterized by a density profile that is assumed to 
take the Navarro-Frenk- White (NFW) form [19]: 

Pnfw = jz f—j n 2 ; (19) 
r{l + r/r s ) 2 
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Fig. 3. Density profiles for four different source models — (a): isotropic model; (b): 
dipole enhancement model (a = 1.0); (c): seven source model (F s = 0.28); (d): dark 
matter halo model (r s = 10 kpc). All four figures are shown in a Hammer- Aitoff 
projection of equatorial coordinates (right ascension right to left). The highest den- 
sity in each panel corresponds to the lightest (red) regions, the lowest density to 
the darkest (blue) regions. 



where p Q is a dark matter density parameter, r is the distance from the center 
of the halo, and r s is a critical radius. For our source model, we will consider 
the contribution of only the four closest significant dark matter halos: the 
Milky Way, M31, LMC, and M33. We will assume that p Q is the same for 
all four sources and that r s scales with the cube root of the luminosity, Ls. 
Thus the Milky Way will have: r S Mw — 10.0 kpc, LMC will have r S LMC — 
0.3 • r S] Mw = 3.0 kpc, M31 will have r Si M3i = 1-5 • t Si mw = 15-0 kpc, and M33 
will have r SjM 33 = 0.4 • r SjM w = 4.0 kpc. 

We now calculate the information dimension, Dj, for each of the our four mod- 
els. Since these are smooth distributions with no statistical fluctuations, we 
will only use one, computationally limited value for the number of declina- 
tional bands for each model of N$ = 1080 (i.e. A9 = I ), which implies: 



(20) 



Using equations 14 and 20 we can now calculate D\ for each of the four models. 
The results are in table 2 column 1. Note that three of the four values of D\ 
exceed the analytical limit of 2 for a 2-D surface. 
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Intrinsic D\ values for the different source models that we examined in 

figures 3 and 5. 
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Isotropic 


2.035 


1.967 


Dipole Enhancement 


2.007 


1.945 


Seven Source 


2.033 


1.946 


Dark Matter Halo Model 


1.999 


1.978 



Table 2 

The estimated values of D\ with N$ = 1800 for the four proposed source models for 
the entire sky independent of any real detector's aperture (1) and for the four source 
models superimposed on the estimated aperture of an air-fluorescence detector at 
40° N. These values are mathematical descriptors of a data set that are whose 
number of significant digits are determined by how extensively each bin is sampled 
(in this case four digits). In the case of a real data set with a finite number of 
the observations, the number of significant digits is constrained by the fluctuations 
inherent to the sample size 



If we consider the analytic limit for the isotropic case, we have: 
(NsfPi log P ; 



Di = -4- 



71 log N s 



(21) 



If we then substitute in equation 20 we get: 
21ogiV ( 5 + log4/7r 



; lim DANs) = 2. 

logiV^ 'AT^oo U dJ 



(22) 



The reason the we obtain values greater than 2 is because we are working with 
a finite number of elements on a surface where the total area does not equal 



4-2 Exposure-Dependent Source Descriptions 



For the purpose of this paper, we assume a hypothetical air-fluorescence de- 
tector located at 40° N. This analysis assumes an isotropic distribution for 
the azimuthal component of the arrival directions and a zenith angle distri- 
bution that remains constant in time. This is what one would expect for a 
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(a) Zenith Angle (Degrees) (b) Sidereal Time (Hours) 

Fig. 4. (a): The distribution of zenith angles for a monocular air- fluorescence de- 
tector; (b): the distribution sidereal times for a air- fluorescence detector located at 
40° N in a desert locale (right) 

detector with 360° coverage with identical detector units and stable atmo- 
spheric conditions. The acceptance of our detector can thus be defined by two 
distributions: zenith angle and sidereal time which are shown in figure 4. The 
zenith angle distribution is characterized by 100% acceptance until ~ 50°, at 
which point it drops off dramatically due to the lack of a well-defined pro- 
file to assist monocular reconstruction. The sidereal time distribution is the 
combination of the seasonal availability of dark, moonless sky at 40° N and 
seasonal climatic changes in a desert locale (the rainy season is assumed to 
extend from February to May which results in a loss of ~ 30% of exposure). 

By defining acceptance this way, we can calculate the exposure of the detec- 
tor. By also considering the finite, asymmetric angular resolution, we can then 
superimpose the detector exposure upon the various source models that we 
previously examined (figure 3 and table 2 column 1). We then obtain an ef- 
fective detector response for the air fluorescence detector for each of our four 
source models. The results are shown in figure 5. It should be emphasized 
that effective detector effect is due to the combined effect of asymmetric sky 
coverage and angular resolution smearing. A remarkable consequence of this 
combination is the possibility that point sources can take on an apparently 
asymmetric shapes due to preferential orientations of the plane of reconstruc- 
tion for events arriving from a specific location in the physical sky. 

We can now determine the value of D\ using the same method as before. The 
results are shown in table 2 column 2. It is interesting to note that the dark 
matter halo source model now has a larger value for Di than the isotropic 
source model. By looking at figure 5, one can verify that the superposition 
of the detector exposure and source models actually yields a more uniform 
apparent distribution for the dark matter halo source model than it does for 
the isotropic source model. 
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Fig. 5. Effective detector response for the four different source models — (a): 
isotropic model; (b): dipole enhancement model (a = 1.0); (c): seven source model 
(F s = 0.28); (d): dark matter halo model (r s = 10 kpc). All four figures are shown 
in a Hammer-Aitoff projection of equatorial coordinates (right ascension right to 
left). The highest density in each panel corresponds to the lightest (red) regions, 
the lowest density to the darkest (blue) regions. 

5 Calculating Di for Finite Event Samples 

So far, we have only considered calculating D\ for smooth distributions. From 
an experimental standpoint, it is very difficult to collect enough data to obtain 
a smooth distribution. This is especially true for UHECRs. In order to make a 
measurement of Di, we must first determine what value(s) we should assign to 
Ng. A reasonable approach is to assign a scale length to our sample. Choosing 
A9 = 0.5°, which approximately reflects the lowest value that can be obtained 
from <t 2 in equation 6, yields N$ = 360. 

However, it can be beneficial to actually calculate D\ for a range of values 
around N$. In figure 6a, we display the values Di(N$) over the range N$ = 
[354, 360] for two separate simulated sets where Ashower = 500. These sets 
yield very similar values for D\ over the full range of values for N s . However, 
if we examine figure 6b we can see that for the same two finite samples, 
the fractional difference between values of D\ can fluctuate substantially even 
over small intervals of A^. While these fluctuations are typically much smaller 
than the difference in D\ values between any two sets, we will take D\ to be 
<Di(N s )> for the interval: Ns = [354, 360] in order to minimize the statistical 
fluctuations in D : between individual sets. This range was chosen to optimize 
our computational ability. 

To account for the fact that we no longer have a smooth distribution for the 
calculation of the values of P ; , we refer back to equation 13 to see how to 
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Fig. 6. (a): D\ over a range of values of N$ for two similar finite event sets that 
have similar values for Dj. (the dots indicate the values D\ for the first set while 
the stars indicate the values D\ for the second set); (b): the fractional difference in 
D\ for the same two sets over the range of Ng values. 

calculate D\ for a finite set of elements. This only requires us to determine a 
value for iVrjist- We can set this value based upon what value we wish for 
in combination with equation 12 (i.e. <rii>= ^ A ^ 2 where An is the fractional 
Gaussian fluctuation of a bin with n ; =<n ; >). Then, 



Dist 



N s 



hower 



(23) 



If we then combine equations 9, 13, and 23; we obtain: 

r 3 



Pi(N 5 



Hi 



71" 



(24) 



We can then calculate D\ from equation 14: 

/ 1 N \ 

Dl = -r^E Fi W 1 °g ?i ( iV 4 N * = [3H360]. (25) 
\ log iV 5 i=1 I 

Thus, If we want <ri\>= 500 and iVshower = 500, we find that ./Voist — 1-65 x 
10 5 . If A^ Showcr = 2000, we have iV Dis t ~ 4.1 x 10 4 . 

We now consider two cases: finite event sets with 500 events and finite event 
sets with 2000 events. These sets will have the angular resolution character- 
istics described in equations 5 and 6. The exposure will be modeled via the 
zenith angle and sidereal time distributions shown in figure 4. Figure 7 con- 
tains examples of event sets with all four source models and iV Shower = 500 
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Fig. 7. Simulated 500 event distributions for four different source models — (a): 
Isotropic model; (b): dipole enhancement model (a = 1.0); (c): seven source model 
(F s = 0.28); (d): dark matter halo model (r s = 10 kpc). All four figures are shown 
in a Hammer-Aitoff projection of equatorial coordinates (right ascension right to 
left). The highest density in each panel corresponds to the lightest (red) regions, 
the lowest density to the darkest (blue) regions. 




Fig. 8. Simulated 2000 event distributions for four different source models — (a): 
Isotropic model; (b): dipole enhancement model (a = 1.0); (c): seven source model 
(F s = 0.28); (d): dark matter halo model (r s = 10 kpc). All four figures are shown 
in a Hammer-Aitoff projection of equatorial coordinates (right ascension right to 
left). The highest density in each panel corresponds to the lightest (red) regions, 
the lowest density to the darkest (blue) regions. 

and figure 8 contains examples of events sets with all four previously described 
source models and A^ower = 2000. 

In figures 7 and 8 one can see the that these distributions of arrival directions 
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( c ) D,-lsotropic vs. Dark Matter Halo 

Fig. 9. Comparison of distributions of D\ values for ~ 500 sets of 500 events between 
the isotopic source model (shaded) and the other three source models — (a): Dipole 
enhancement model vs. isotropic model (a = 1.0); (b): Seven source model vs. 
isotropic model (F s = 0.28); (c): Isotropic model vs. Dark matter halo source model 
(r s = 10 kpc). The vertical line corresponds to the value of D\ for the "real" 500 
event sample. In all cases the distributions of D\ values fit well to a Gaussian curve. 



have a far greater degree of statistical fluctuation than the smooth distribu- 
tions shown in figure 5. Because of the fluctuations in our simulated event 
samples, the value of D\ varies significantly (see figures 9 and 10) from one 
simulated set to the next. In figures 9 and 10, we examine the distribution of 
D\ values for ~ 500 sets (of 500 and 2000 events respectively). We see that 
the distribution of 500 event samples have both a lower mean value and larger 
width than the 2000 event samples. 
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( c ) D,-lsotropic vs. Dark Matter Halo 



Fig. 10. Comparison of distributions of D\ values for ~ 500 sets of 2000 events 
between the isotopic source model (shaded) and the other three source models - 
(a): Dipole enhancement model vs. Isotropic model (a = 1.0); (b): Seven source 
model vs. Isotropic model (F s = 0.28); (c): Isotropic model vs. Dark matter halo 
source model (r s = 10 kpc). The vertical line corresponds to the value of D\ for 
the "real" 2000 event sample. In all cases the distributions of Dj values fit well to 
a Gaussian curve. 

6 Application to Anisotropy Analysis 

We now need to develop a scheme by which we can apply fractal dimensionality 
analysis to a real data set. In the case of a real data set, we will be dealing with 
only one value of D\. By itself D\ is insufficient to characterize the data set; 
D\ fluctuates a great deal due to variation in iVshower- However, a comparison 
between the value of D\ for a real event sample and a distribution of D\ 
values (with the same N ShowCT and iV Dist values as the real data) for a series of 
simulated data sets of a given source model does provide a viable measurement 
of anisotropy. 

We can demonstrate this by considering a single simulated event sample gen- 
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erated with an isotropic source model. We will suppose this sample to be our 
"real" data. We consider the isotropic simulated event samples shown in fig- 
ures 7a and 8a. We will once again stipulate that <rii>= 500 which means 
that in the case of N S hower = 500, we have: iV Dist ~ 1.6 x 10 5 which leads to 
Di = 1.89715 and in the case of N S hower = 2000, we have: N Dist ~ 4 x 10 4 
which leads to Di = 1.93920. 

In figures 9 and 10, we demonstrated that for a fixed scaling parameter the 
distribution of D\ values for a large number of simulated sets fits well to a 
Gaussian curve. By establishing the relationship between the distribution of 
Di values and the scaling parameter in each model, we can establish a 90% 
confidence interval on the scaling parameter for that model. 

6. 1 Dipole Enhancement Source Model 

In the case of the dipole enhancement source model in equation 17, the scaling 
parameter is a. By varying a between —1 and 1, we develop a curve which will 
show the relationship between Di and a. By considering the actual value of 
Di for the "real" data set, we then establish a nominal value for a and a 90% 
confidence interval. The results for both N S hower — 500 and N Shower = 2000 
are shown in figure 11. In the case of our simulated isotropic set with Nshower = 
500, a = 0.02 ± 0.21 with a 90% confidence interval of [-0.29,0.36]. In the 
case of our simulated isotropic set with Nshower — 2000, a = 0.075 ± 0.085 
with a 90% confidence interval of [—0.065, 0.24]. These results tell us that our 
"real" data can, at the la level of up to 0.2. We notice that in figure 11, 
there are potentially two suitable intervals of a that possess similar values of 
D\. This emphasizes the necessity of making density plots to assure by visual 
inspection that the appropriate interval is chosen. 

6.2 Seven Source Model 

In the case of the seven source model, the scaling parameter is F s , the fraction 
of the total event sample that is produced by the discrete sources. By vary- 
ing F s between and 0.40, we develop a curve which shows the relationship 
between D\ and F s . By considering the actual value of D\ for the "real" data 
set, we then establish a nominal value for F s and a 90% confidence upper limit 
for F s . The results for both N shower — 500 and N shower — 2000 are shown in 
figure 12. In the case of our simulated isotropic set with N Shower = 500, we 
have a 90% confidence upper limit of F s = 0.16. In the case of our simulated 
isotropic set with N shower = 2000, we have a 90% confidence upper limit of 
F s = 0.04. This tells us that our "real" data can have, at most, 16% (for 500 
events) or 4% (for 2000 events) of these events coming from the seven sources. 
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(c) a (d) 

Fig. 11. Study of the dependence of D\ upon a for a dipole enhancement source 
model— (a): N Shower = 500 a = .02 ± .21; (b): N Shower . = 500 (zoomed); 90% 
confidence interval: [-0.29,0.36]; (c): N Shower = 2000; a = .075 ± .085; (d): 
N shower = 2000 (zoomed); 90% confidence interval: [—0.065,0.24]. In each case, 
the solid solid horizontal line indicates the value of D\ for the simulated isotropic 
data set. The vertical lines in (b) and (d) indicate the projection of the nominal 
value and 90% and 95% confidence intervals of a on the x-axis. 



6.3 Dark Matter Halo Source Model 



In the case of the dark matter halo source model in equation 19, the variable 
parameter is r s , the critical radius in the NFW profile [19]. By varying r s 
between 5.0 kpc and 15.0 kpc we develop the curve which will demonstrate 
the dependence of D\ upon r s . We can then show that the dark matter halo 
source model can be rejected with at level > 3.6<r for Nghower — 500 and at 
a level of > 7.0a for N Shower = 2000 for the full range of hypothesized values 
for r s . The results are shown in figure 13. 
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Fig. 12. Study of the dependence of D\ upon F s for the seven source model — (a): 
N Sh ower — 500; (b): TV shower — 500 (zoomed); 90% confidence upper limit— 0.16; (c): 
N shower = 2000; (d): N shower = 2000 (zoomed); 90% confidence upper limite=0.04;. 
In each case, the solid solid horizontal line indicates the value of D\ for the simulated 
isotropic data set. The vertical lines in (b) and (d) indicate the projection of 90% 
and 95% confidence upper limits of r s on the x-axis. 



7 Discussion 



Fractal dimensionality has several advantages over conventional anisotropy 
techniques. First of all, it naturally accommodates angular resolution. This 
is extremely important when considering event sets with asymmetric errors, 
when analysing event sets with variable values for the angular resolution (e.g. 
dependent on energy or geometry), or when combining multiple data sets from 
different detectors for a single analysis. 

Another advantage that fractal dimensionality possesses is the ability to ac- 
commodate any aperture. Because this method makes a relative comparison 
between a sample and simulations using the same aperture, the physical aper- 
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Fig. 13. Comparison of the distribution of D\ values for a dark matter halo source 
model with a full range of hypothesized values for r s — (a): N 'shower = 500, the dark 
matter halo model can be rejected at a level > 3.6a; (b): Nghower = 2000, the 
dark matter halo model can be rejected at a level > 7.0a. In each case, the solid 
horizontal line indicates the value of D\ for the simulated isotropic data set. 

ture is simply folded into the analysis. This once again allows the combination 
of data from multiple detectors with very different apertures. It also allows the 
analysis of extremely complicated apertures without the need to include nor- 
malizing factors that needlessly complicate the predictability of the Poisson 
fluctuations in the data sample. 

Perhaps the most striking feature of the fractal dimensionality method is that 
it only requires a single measurement of one's data. While fractal dimen- 
sionality will not always provide better statistical significance than a direct 
measurement for a specific anisotropy, the fact that one considers only a sin- 
gle measurement of the data, for any number of potential anisotropic models, 
provides one the means to simultaneously reject or accept all of those models 
without the ensuing statistical penalty. 

A possible way of increasing the sensitivity of the fractal dimensionality method 
is by considering the general case of case of D q . By varying the value of q to 
something other than 1, it might be possible to increase the sensitivity of this 
method to various anisotropics. 

The fractal dimensionality method does have some weaknesses. First and fore- 
most is the potential for multiple solutions as was demonstrated in the dipole 
source model above. The method cannot be applied blindly. It requires a care- 
ful inspection of both the data sample and the simulated samples in order to 
resolve possible ambiguities. Another drawback is the amount of computation 
required to calculate D\ for a large number of simulated data sets. Producing 
just the plots in figure 11c consumed the equivalent of ~ 1000 CPU hours on 
1GHz Athlon machine. Another limitation is the potential for different source 
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models to effectively cancel each other out and yield a potentially deceiving 
value of Di that resembles that of an isotropic source. One solution for this 
is to separately consider the value of D\ for different celestial regions in the 
data. Of course, this will incur a statistical penalty. 

Two particularly useful roles for fractal dimensionality analysis are (a): as 
a first test to ascertain if a sample possesses the same heterogeneity as the 
expected isotropic background and (b): as a wholly independent observable 
that can be used to confirm the results of a direct measurement. 

Today the HiRes detector continues to accumulate data. Soon the Auger de- 
tector will be acquiring data at seven times the rate of HiRes. As the number 
of observed UHECRs continues to rise, fractal dimensionality will grow ever 
more effective in its ability to discern the between potential source models in 
the continuing effort to solve the mystery of UHECRs. 
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